STATISTICAL STABILITY IN TIME REVERSAL 
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Abstract. When a signal is emitted from a source, recorded by an array of transducers, time reversed 
and re-emitted into the medium, it will refocus approximately on the source location. We analyze the refo- 
cusing resolution in a high frequency, remote sensing regime, and show that, because of multiple scattering, 
in an inhomogeneous or random medium it can improve beyond the diffraction limit. We also show that 
the back-propagated signal from a spatially localized narrow-band source is self-averaging, or statistically 
stable, and relate this to the self-averaging properties of functionals of the Wigner distribution in phase 
space. Time reversal from spatially distributed sources is self-averaging only for broad-band signals. The 
array of transducers operates in a remote-sensing regime so we analyze time reversal with the parabolic or 
paraxial wave equation. 
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1. Introduction. In time reversal experiments a signal emitted by a localized source 
is recorded by an array and then re-emitted into the medium time-reversed, that is, the tail 
of the recorded signal is sent back first. In the absence of absorption the re-emitted signal 
propagates back toward the source and focuses approximately on it. This phenomenon 
has numerous applications in medicine, underwater acoustics and elsewhere and has been 
extensively studied in the literature, both from the experimental and theoretical points of 
view fl|, H, P, ||, |o], || ||, H. Recently time reversal has been also the subject 
of active mathematical research in the context of wave propagation and imaging in random 
media [^0, ]?], |[ |], ^2) . A schematic description of a time reversal experiment is presented 
in Figure 



XL 




Fig. 1.1. A pulse propagates toward a time reversal array of size a. The propagation distance L 
is large compared to a. The ambient medium has a randomly varying index of refraction with a typical 
correlation length that is small compared to a. The signal is time reversed at the array and sent back into 
the medium. The backjaropagated signal refocuses with spot size \L/a e , where a e is the effective aperture 
of the array (Section 3.1 ). 
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For a point source in a homogeneous medium, the size of the refocused spot is approx- 
imately XL /a, where A is the central wavelength of the emitted signal, L is the distance 
between the source and the transducer array and a is the size of the array. We assume here 
that the array is operating in the remote-sensing regime a <ti L. Multiple scattering in a 
randomly inhomogeneous medium creates multipathing, which means that the transducer 
array can capture waves that were initially moving away from it but get scattered onto it 
by the inhomogeneities. As a result, the array captures a wider aperture of rays emanating 
from the original source and appears to be larger than its physical size. Therefore, somewhat 
contrary to intuition, the inhomogeneities of the medium do not destroy the refocusing but 
enhance its resolution. The refocused spot is now XL/a e , where a e > a is the effective 
size of the array in the randomly scattering medium, and depends on L. The enhancement 
of refocusing resolution by multipathing is called super-resolution 0. The time reversed 
pulse is also self averaging and refocusing near the source is therefore statistically stable, 
which means that it does not depend on the particular realization of the random medium. 
There is some loss of energy in the refocused signal because of scattering away from the 
array but this can be overcome by amplification, up to a point. 

The purpose of this paper is to explore in detail the mathematical basis of pulse sta- 
bilization, beyond what was done in Q. We want to explore in particular in what regime 
of parameters statistical stability is observed in time reversal. We show here that for high 
frequency waves in a remote sensing regime, spatially localized sources lead to statistically 
stable super-resolution in time reversal, even for narrow-band signals. We also show that 
when the source is spatially distributed, only for broad-band signals do we have statistical 
stability in time reversal. The regime where our analysis holds is a high frequency one, more 
appropriate to optical or infrared time reversal than to ultrasound, sonar or microwave radar. 
In this regime we can make precise what spatially localized or distributed means (see Section 



3.1). The numerical simulations in J7| and |§J], which are set in ultrasound or underwater 
sound regime, indicate that time reversal is not statistically stable for narrow-band signals 
even for localized sources. Only for broad-band signals is time reversal statistically stable 
in the regime of ultrasound experiments or sonar. 

If the aperture of the transducer array is small a/ 'L -C 1, the Fresnel number L/(ka 2 ) 
is of order one, and the random inhomogeneities are weak, which is often the case, we may 
analyze wave propagation in the paraxial or parabolic approximation . The wave field is 
then given approximately by 

(1.1) u(t, x, z) = i J e^'^-^iz, x; Lo/c )duo 

where the complex amplitude ip satisfies the parabolic or Schrodinger equation 

(1.2) 2ihp z + A x iP + k 2 (n 2 - l)tp = 0. 

Here x = (x,y) are the coordinates transverse to the direction of propagation z, the wave 
number k = uj/cq and n(x, z) — cq/c(x,z) is the random index of refraction relative to a 
reference speed cq. The fluctuations of the refraction index 

(1-3) <r/,(y,|) = n 2 (x,z)-l 

are assumed to be a stationary random field with mean zero, variance a 2 , correlation length 
I and normalized covariance with dimensionless arguments 



(1.4) 



i?(x, z) = E{fi(x + x', z + z')M x 'i z ')}- 
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A convenient tool for the analysis of wave propagation in a random medium is the 
Wigner distribution |l9|, p8| defined by 



(1.5) 



W{z, x, p) = — / e^(x - L)^(x + ^,z)dy 



where d = 1 or 2 is the transverse dimension and the bar denotes complex conjugate. The 
Wigner distribution may be interpreted as phase space wave energy and is particularly well 
suited for high frequency asymptotics and random media p8| . The quantity of principal 
interest in time reversal, the time-reversed and back-propagated wave field, can be also 



expressed in terms of the Wigner distribution (see Section |3.l| ) . The self-averaging properties 
of the back-propagated field are related to the self-averaging properties of functionals of the 
Wigner distribution in the form of integrals of W over the wave numbers p. 

In the next Section we introduce a precise scaling that corresponds to (a) high frequency, 
(b) long propagation distance, (c) narrow beam propagation, and (d) weak random fluctua- 
tions. In the asymptotic limit where the small parameters go to zero the Wigner distribution 
satisfies a stochastic partial differential equation (SPDE), a Liouville-Ito equation, that has 
the form 

(1.6) dW(z, x, p; A;) = • V X W + ^-A p W^j dz - ^V p W • dB(x, z) 
where B(x, z) is a vector- valued Brownian field with covariance 

(1.7) E{B i ( Xl ,z 1 )B j (^ 2 ,z 2 )} = - ( ^^g"^^ ) Zl A Z2 ' 
where Z\ A z 2 = min{zi, z 2 }, and in the isotropic case 

(1.8) D = -^°p, i? ( x ) = r R{*,s)ds. 



In Section 2.5 we analyze this SPDE in the asymptotic limit of small correlation length 
for B(x, z) in the transverse variables x, and show that W(z, x, p; k)'s with different wave 
vectors p are uncorrelated. From this decorrelation property we deduce that for localized 
sources the time-reversed, back-propagated field is self-averaging, even for narrow-band sig- 
nals. For distributed sources it is self-averaging only for broad-band signals. We show in 
detail in Section ^ how the asymptotic theory is used in time reversal. In Appendix |A| we 
introduce other scalings which lead to the same averaged SPDE but we do not analyze them 
in detail. 

Throughout the paper we define the Fourier transform by 

/(k) = ydxe- ikx /(x) 

so that 
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Sloan Fellowship and ONR grant N00014-02- 1-0089. K. Solna by NSF grant DMS-0093992 
and ONR grant N00014-02-1-0090. 
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2. Scaling and asymptotics. 

2.1. The rescaled probl em . To carry out the asymptotic analysis we begin by rewrit- 
ing the Schrodinger equation (1.2) in dimensionless form. Let L z and L x be characteristic 
length scales in the propagation direction, as, for example, the distance L between the source 
and the transducer array for L z and the array size a for L x . We introduce a dimensionless 
wave number k' — k/ko with ko = ^o/co and luq a central frequency. We rescale x and z by 



x = i x x', z = L z z' and rewrite (1.2) in the new coordinates dropping primes: 



(2-1) 



2ik 



dtp 
dz 



knL': 



-Atp + k koL z a/i 



I 



zL z 
I 



ip = 0. 



The physical parameters that characterize the propagation problem are: (a) the central 
wave number ko, (b) the strength of the fluctuations <r, and (c) the correlation length I. We 
introduce now three dimensionless variables 



(2-2) 



S-— E — — - — 

L x L z ' fco^ 



which are the reciprocals of the transverse scale relative to correlation length, the recip- 
rocal of the propagation distance relative to correlation length, and the central wave 
length relative to the correlation length. We will assume that the dimensionless parameters 
7, <7, e and 5 are small 



(2-3) 



7<Cl; cr < 1; <5 < 1; e < 1 



This is a regime of parameters where super-resolution phenomena can be observed. 
To make the scaling more precise we introduce the Fresnel number 



(2.4) 



L z P 
k L^ ^ e 



We can then rewrite the Schrodinger equation (2.1) in the form 



(2.5) 



k^S x z 

£ L ' Z £ 



provided that we relate e to a and 5 by 

(2.6) £ = a 2 'H 2 '\ 



One way that the asymptotic regime (2J3) can be realized is with the ordering 



(2.7) 



< £ < 5 < 1 



and 7 <§; cr 4 / 3 <5~ 2 / 3 , corresponding to the high-frequency limit. We see from the scaled 
Schrodinger equation (2.5) that this regime can be given the following interpretation. We 
have first a high frequency limit #—►(), then a white noise limit e — > 0, and then 
a broad beam limit 5 — ► 0. We will analyz e in detail and interpret these limits in the 
following Sections. Another scaling in which (2.3) is realized ise <C 6 <C$ « 1. This 
is a regime in which the white noise limit is carried out first, then the high frequency limit 
and then the broad beam limit. We do not analyze this case here. Additional comments on 
scaling are provided in Appendix 
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It is instructive to express the constraints ( (2.6| ) and (2.7) in terms of the dimensional 
parameters of the problem. First, both the size of the transverse scale L x and the propagation 
distance L z should be much larger than the correlation length I of the medium. Moreover, 
( ^ ) implies that the longitudinal and transverse scales should be related by 

> 1 



L z ( * N 1/3 



so that we are indeed in the beam approximation. The first inequality in (2.7) implies that 



« VM = -=, 

L X y/j 

and with the above choice of L z this implies that 



7' 



3 / 2 1 



< — < 



a 2 I a 2 ' 

2.2. The high frequency limit. A convenient tool for the study of the high frequency 
limit, especially in random media, is the Wigner distribution. It is often used in the context 
of energy propagation [ fl9[ p8| but it is also useful in analyzing time reversal phenomena 
(2[ H Let 4>e(x) be a family of functions oscillating on a small scale 9. The Wigner 
distribution is a function of the physical space coordinate x and wave vector p defined as 



(2.8) Wo(x,p) = |^e^ e (x-^)^(x+^). 

The family Wg is bounded in the space of Schwartz distributions 6>'(R d xl°) if the functions 
4>g are uniformly bounded in L 2 (H d ). Therefore there exists a subsequence 9 k — > such that 
W$ k converges weakly as k — > oo to a limit measure W^(x, p). This limit W(x, p) is non- 
negative and is customarily interpreted as the limit phase space energy density because 

(2.9) |0 Sfc (x)| 2 ^ yV(x,p)dp as 9^0 

in the weak sense. This allows one to think of W(x, p) as a local energy density. 

L et W g(z, x, p) be the Wigner distribution of the solution ip of the Schrodinger equa- 



tion (2.5), in the transversal space-variable x. A straightforward calculation shows that 



Wg(z,x,p) satisfies in a weak sense the linear evolution equation 



dWa p 



kS 



^(P-H)-^(P + If) rfq 

In the limit 9 — > the solution converges weakly in S', for each realization, to the (weak) 
solution of the random Liouville equation 

The initial condition at z = is VT(0, x, p) = Wi (x, p), the limit Wigner distribution of the 
initial wave function. 
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2.3. The white noise li mit. In this Section we take the white noise limit e — > in 
the random Liouville equation ( [2.1l| ) whose solution we now denote by W e . We can do this 
using the asymptotic theory of stochastic differential equations and flows ^, |T[ pp| a s 
follows. Using the method of characteristics, the solution of the Liouville equation ( 2.11 ) 
may be written in the form 

W £ (t, x, p) = Wj (X E (i; x, p) , P £ (t; x, p)), 

where the processes X e (t; x, p) and P e (i; x, p) are solutions of the characteristic equations 

dX F 1 dP e k „ fX e z 



dz k 1 dz 2y/e 

with the initial conditions X e (0) = x and P e (0) = p. The asymptotic theory of random 
differential equations with rapidly oscillating coefficients implies that, under suitable condi- 
tions on /x, in the limit e — > the processes X e , P £ converge weakly (in the probabilistic 
sense), and uniformly on compact sets in x, p to the limit processes X(t), P(t) that satisfy 
a system of stochastic differential equations 

dP = --dB(z), dX = ~Pdz, X(0) = x, P(0) = p. 
2 k 

The random process B(z) is a Brownian motion with the covariance function 
(2.12) E{B i (z 1 )B j (z 2 )}= - d Q^ ds Zl A z 2 



= <% (--Ro(O)) «i A 

in the isotropic case, where 



•-2- 



(2.13) Ro(x) = / R(x,s)ds 

is a function of |x|. This implies that the average Wigner distribution We {z, x, p) = 
E {W e {z, x, p)} converges as e — > uniformly on compact sets to the solution of the 
advection-diffusion equation in phase space 

(2.14, ^ + £. V ^. = ^A pW '<.) 

oz k 2 

with the initial data (0, x, p) = W/(x, p). Here the diffusion coefficient D is given by 

(2.15) D = -^. 

4 

The one-point moments E { [W e (^, x, p)]^} converge as e — > to the functions (z, x, p) 
that satisfy the same equation ( |2.14[) but with the initial data (0, x, p) = [Wo(x, p)] N ■ 
This is similar to the spot dancing phenomenon where all one-point moments are gov- 
erned by the same Brownian motion. In particular we have that 

W^(z,x,p)^ [^«(z,x,p) n2 
so that the process W £ does not converge to a deterministic one, in the strong sense pointwise. 
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2.4. Multi-point moment equations. As in the previous Section we may also study 
the white noise limit e — > of the higher moments of W e (z, x, p) at different points 

^W(z,x\...,x w ,p\...,p w )=^{[^(z,x\p 1 )P-...-[^(z,x w ,p^)r}. 

Here the points (x m ,p m ) are all distinct, (x",p n ) ^ (x m ,p m ). We may account for mo- 
ments that have different powers of W e at different points by taking different powers Tj of 
W ix-'.p':. 

We now consider the joint process (X e (z; x m , p m ), P e (^; x m , p" 1 )), m = 1, ...,N. As 
e — ► it converges to the solution of the system of stochastic differential equations 

u N d / Y' n — X" \ 1 
(2.16) dPT = - J E E^ 1 dB i {z) ' = ~k Pmdz ' 

with the initial conditions 

X m (0)=x m , P m (0)=p m . 
The d-dimensional Brownian motions B m , m = 1, . . . , N have the standard covariance tensor 

E{B™( Zl )B?(z 2 )} =5 mn 6 i jZ 1 Az 2 , i, J = 1, . . . , d, m, n = 1, . . . , N. 
The symmetric tensor er^ (x) is determined from 

JV 



We assume that equation ( |2.17| ) has a solution that is differentiable in x, which is compatible 
with the fact that the matrix on the right is, by Bochner's theorem, non-negative definite. 
The moments converge as e — > to the solution of the advection-diffusion equation 

(2ig) ^^t^.w^^ 

m=l m—1 

_ — V V d 2 i?o((x"-x m )/<5) d 2 W^ 

71, 771=1 l 7 J — l J J 

n>m 

with the initial data 

W^(0, xi, , . . . , x w , p 1 , . . . , p w ) = [Wr(x\ p 1 )]^ [W/(x w , p w )F ■ 
From ( |2.18| ) we can calculate moments of functionals of W e of the form 

W e ,${z) = J W e (z,x, p)</>(x,p)dxdp. 

For example, as £ — ► we have that 

E {[W e ^{z)] 2 } -> y M /(2) (z,xi,pi,x 2 ,p 2 )0(xi,pi) ( Tj(x2,P2)dxidpidx2dp2. 

A convenient way to deal with not only the limit of Appoint moments but with the 
full limit process W(z, x, p), at all points x, p simultaneously, is provided by the theory of 
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stochastic flows [£3|. For this we need to show that W e {z, x, p) converges weakly (in the 
probabilistic sense) as e — > to the process W(z,x, p) that satisfies the stochastic partial 
differential equation 



(2.19) 



dW s 



fc 2 .D 



5 



dz - -V p Ws ■ dB(- 



Here the Gaussian random field B(x, z) has the covariance 



J B{B i (xi,zi)B i (x 2 ,Z2)} 



a 2 fi ((xi-x 2 )) 

dxidx-j 



Z\ A z 2 . 



We call equation (2.19) the Liouville-Ito equation. It allows us to treat all equations of the 
form (2.18) simultaneously and is a convenient tool for simulation and analysis. The dimen- 
sionless wave number k can be scaled out of (2.19) by writing W(z,x, p; k) = W(z, x, §; 1) 
so that we need only consider (2.19) with k = 1. We will use this scaling in Section [O. 

Note that unlike the single Brownian motion (2.12) that governs the evolution of one- 
point moments, the Brownian field that enters the SPDE (2.19) depends explicitly on the 
dimensionless correlation length S in the transverse direction. Therefore the limit process 
also depends on 6 and we denote it by Wg. 

2.5. Statistical stability in the broad beam limit. We will now consider the limit 
8 — > of the process Wg(z,x,p) when the transverse dimension d > 2. We are particularly 
interested in the behavior of functionals of Wg as 6 — > 0. The analysis of one-point moments 



in Section 2.2 showed that they do not depend on 5 and are governed by a standard Brownian 
motion. Therefore the process Wg does not have a pointwise deterministic limit. However, 
we will show that functionals of Wg become deterministic in the limit S — > 0. We refer to this 
phenomenon as statistical stabilization and give conditions for it to happen. Stabilization 
plays an important role in time reversal, imaging and other applications, as discussed in the 
Introduction. 

Theorem 2.1. Assume that<p(p) is a smooth test function of rapid decay, the transverse 
correlation function i?o( x ) has compact support, the initial Wigner distribution W/(x, p) is 
uniformly bounded and Lipschitz continuous, and the transverse dimension d > 2. Define 



(2.20) 

Then 
(2.21) 



Ig^{z,^) 



Wg(z,x,p)(j)(p)dp. 



lim£{/! (2,x)} =£ 2 {/^(z,x)} 



where E {/^(^x)} is independent of 5. 

The assumption of compact support for i?o( x ) is not essential but simplifies the proof. 
We have already noted that the Wigner distribution Wg itself does not stabilize. However, 



(2.21) implies that 

(2.22) lim Var = lim E - E 2 = 0. 

Therefore, any smooth functional of the form ( 2.20| ) stabilizes in the limit 5 — > 0, that is, 

(2.23) J 5)0 w E{Ig i<t> }, 
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in mean square, and the expectation of Ig^ does not depend on 5. We prove Theorem 2.1 
in Appendix [b|. 

In the applications of the asymptotic theory to time reversal we need not only functionals 
Is. 4, of the form ( 2.20| ) but also of the form 



(2.24) 



Js(z,x) = J Ws(z,x, p)dp. 



We need to show that such functionals are well defined with probability one and to analyze 
their behavior as 6 — > 0. This is done in the following theorem. 

Theorem 2.2. Under the same hypotheses of Theorem 2.1 and with a non-negative 
initial Wigner distribution Wj > 0, the functional J$ is bounded, continuous and non- 
negative with probability one. In the limit 5 — > we have 



(2.25) 



Jim E { Jj (z, x) } =E 2 {J s {z,x)} 



where E {J$(z,ic)} does not depend on 8. 

The proof of this theorem is given in Appendix! 



What is important in both Theorems 2.1 and 2.2 is that we do integrate over the wave 



numbers p because there is no pointwise stabilization. In time reversal applications, as in 
section 3.1, we actually need Theorem 2.2 when the integration is only over a line segment in 



p sp ace, and the dimension of the latter is d > 2. Its proof follows from the one of Theorem 
2.2. 



3. Application to time reversal in a random medium. We will now apply these 
results to the time reversal problem Q described in the Introduction. A wave emitted 
from the plane z — propagates through the random medium and is recorded on the time 
reversal mirror at L. It is then reversed in time and re-emitted into the medium. The back- 
propagated signal refocuses approximately at the source, as shown in Figure 



1.1 



There are 

two striking features of this refocusing in random media. One is that it is statistically stable, 
that is, it does not depend on the particular realization. The other is super-resolution, that 
is, the refocused spot is tighter than in the deterministic case. We discuss these two issues 
in this section. 

3.1. The time-reversed and back-propagated field. We assume that the wave 
source at z = is distributed on a scale a s around a point Xq, that is, 



i/> e (z = 0, x; k) = e J P"-( x - X0 )/%(- 



x 



k). 



where ipo is a rapidly decaying and smooth function of x and k. The width of the source 
a s could be large or small compared to the Fresnel number 9, and this affects the statistical 
stability of the time-reversed, back-propagated field, as we explain in this Section. The 



Green's function, Ge(.z,x;£), solves the parabolic wave equation (2.5) with a point source 
at (x, z) — (£,0). Using its symmetry properties and the fact that time reversal i — ^ — t is 
equivalent to lj — > —to or k — > —k, the back-propagated, time-reversed field on the plane of 
the source has the form 

(3.1) ^?(L,xo,6fc) = 

G e (L, x; x + fc)G,(L, x + ??; 77; fc)e lpo " )/ Vo( — ; -k) X A^)dicdr,. 

0~s 

The complex field amplitude ipg is evaluated at xo + 6£, in the plane z — 0. We scale the 
observation point off Xq by 6 because we expect that the spot size of the refocused signal 
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will be comparable to the lateral spread of the initial wave function. We denote with \A 
the aperture function of the time reversal mirror. It could be its characteristic function, 
occupying the region A in the plane z = L 



Xa(x) 



1, x e A 
0, x 4 A 



or a more general aperture function like a Gaussian. The time reversal mirror is located in 
the plane z = L. 

After changing variables, the back-propagated field is given by 

9-q 



i/>f(L, x , £; k) = 9 d / G e (L, x; x + k)G e {L, x; x + 6 m k)e^M — ^k) XA (^)d^dr, 



= 9 d f G e (L,x + ^,x;fc)G e (L,x + ^,x;A : )e lpoI 'Vo( — ;-fc) Xj4 (x)dxdry. 
J os 

It is now convenient to introduce the Wigner distribution 

r d e ip y 

(3.2) We (z, xo, P',k) = J Ge{z, x - yd/2, x; k)G e (z, x + yO/2, x; k)xA(x)dxdy, 

and express the back-propagated field as 

(3.3) 1>f(L, xo, k) = [ e^-^W e (L, x + ^^ ,p; k)e^M — 5 -k)dpdr,. 

J 2 a s 



The Wigner distribution is scaled here differently from (2.8) because of the way we have 
scaled the source function. 

In the high frequency limit 9^0, W$(z, x, p; k) tends to W(z, x, p; k), which solves the 
random Liouville equation (2.11). Then, in the white noise limit, it solves the Liouville-Ito 
equation (2.19). The mean of W solves (2.14), in the high-frequency and white noise limit, 
with initial data 



(3.4) 

Let 
(3.5) 



W(0,x,p;fc) = 



Xa(x) 
(2ir) d ' 



be the ratio of the width of the source to the Fresnel number and assume that it remains 
fixed as 9 — > 0. In this limit, the time-reversed and back-propagated field is given by 

(3.6) ^ s 0L,x o ,£;fc)= / e J P-(«-")W(L,xo,|)e lpo -^o(r///?;-fc)dpdr, 



e ip <W(L, x , f )[3 d MP(P - Po); ~k)dp. 



Here we have used the scaling W(z,x,p;fc) = W(z,x, 1) in ( 2.1S ) and we have dropped 
the last argument k = 1. 
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3.2. Statistical stability. From the form (|3.6| ) of the back-propagated and time- 
reversed field we see that when (3 = 0(1) (or small), which means that a s is comparable 



to the Fresnel number 9 (or smaller), we can apply the results of Section 2.5 and conclude 
that it is statistically stable or self-averaging in the broad beam limit <5 — > 0. Theorems 2.1 



and 2/2 are exactly what is needed for this. The fact that the initial function (3^4) may 
be discontinuous at the boundary of the set A is not a problem. This is because, we may 
approximate the function \A from above and below by two smooth positive functions, to 
which we may apply Theorems 2.1 and |2.2[ and then use the maximum principle to deduce 
the decorrelation property when the initial data is XA- We have, therefore, 

VA£,x ,£;fc)«(^ s (L,x ,e;fc)} 

in the sense of convergence in probability or in mean square, in the broad beam limit 6 — ► 0, 
for each fixed frequency uj — fcco- Statistical stability of time reversal does not depend on 
having a broad-band signal if the source is localized in space. This is true in the regime of 
parameters reflected by the scaling 8 <C e <C S considered here, which is a high frequency 
regime encountered in optical or infrared applications like ladar. The numerical experiments 
in J7j and || are closer to the regime of ultrasound experiments pl| and in underwater sound 
propagation, which is different from the high frequency regime analyzed here. 



For distributed sources the parameter (3 is large and we cannot apply Theorems 2.1 and 
2.2 to (3.6). It is necessary for statistical stability in this case to have broad-band signals. 
For (3 large the time reversed and back propagated signal in the time domain has the form 



(3.7) 



^ B (L,x ,£,i) 

= (27r) d e l{po ^- kocot) M^/P) I W(L,x 



Po 



k + k 

UJQ + UJ 



)e-* kcot g(-c k) 



c dk 
2tt 



dio 
2^ 



with g(cok) the Fourier transform of the initial pulse relative to the central frequency ojq = 
co&o- This means that we have replaced the actual wave number k by fco + k, or u) by luq + lu, 
with the new oj, the baseband frequency, bounded by f2, the bandwidth, \ui\ < fl < ujq. 
The integration is over the bandwidth [— fi, fl]. This integral is well defi ned with probability 
one and is self-averaging in the broad beam limit 5 — > by Theorem 2.2 and the remark 



following it. We will compute its average in Section 3.4 



3.3. The effective aperture of the array. From the explicit expression for the 



Green's function of (2.14), with k = 1 



C/(z,x,p;x°,p ) 



dwdr 



x exp 



(27T) 

Dz 



2il 



exp (iw ■ (x — x ) + ir ■ (p — p ) — izw ■ p°) 



r + zr ■ w 



and with the time reversal mirror a distance L from the source and xq = 0, it follows from 



(3.6) that 



(3.8) ^ B (z^-k)) 



dpdydw ip< Qd 



(2ir 



,2(/ 



e l P'Wo(/3(p - po); -k) XA (y) exp 



-jw • y — izw ■ — 

k 
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The high-frequency, white-noise limit of the self-averaging time-reversed and back-propagated 
field is therefore given by a convolution 

(3.9) (ip B (L,£,k))=^(;-k)*W(-)(0 
with 

(3.10) W( V ) = W( m L, k) = J^yXA(vk/L) e ~^^\ 
the point spread function, and 

(3.11) -k) = e lp ^Mv/m(-kc ) 

with iPq(j]/(3) the spatial source distribution function and g the Fourier transform of the 
pulse shape function g(t). This notation is consistent with ( |3.7[ ), with the time factor 
e -zk c t om itted, along with the horizontal phase e lkz which cancels in time reversal. We 
have also introduced the refocused spot size with multipathing 



{: '' A ' 2 : a M ~ 77773 " 7771 



I? 

DLk 2 k 2 al 



and the effective aperture a e = a e (L), 



(3,3) 

which we now interpret. 

If the time reversal mirror is the whole plane z = L, then xa = 1 and 

(^ B (L,£;fc)) = V^,-fc). 

In this case the back-propagated field is the source field reversed in time, both in the random 
and in the deterministic case. The point spread function W determines the resolution of 
the refocused signal for a time reversal mirror of finite aperture. Multipathing in a random 



medium gives rise to the Gaussian factor (3.12) whose variance is a M . We can give an 



interpretation of this variance, or spot size, as follows. For a square time reversal mirror of 
size a, the Fourier transform of XA is the sine function so that 

W( m , m ;L,k) = (^)^in(^)sin(^)e-^^)/(^) 

For a deterministic medium (D — 0) the Rayleigh resolution is the distance r)p to the first 
zero of the sine, the first Fresnel zone in either direction, 

2nL XL 

VF = -r— = • 

ka a 



In general, if xa is supported by a region of size a we may define the Fresnel resolution, or 
the Fresnel spot size, by 
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For weak multipathing we have <jm ^> crp and 

W(ri;L,k)~ XA(vk/L) , 

which is the diffractive point spread function whose integral over 77 £ R d is one. If, however, 
we have strong multipathing, <jm "C &f, then we may approximate XA^k/L) by Xa(0) = 
a d in ( 3.10| ), and the point spread function becomes 

W( W X,*)~(^ye-WV(^). 

By writing the variance (spot size) o 2 M in the form ( |3.12 ) we can interpret a e as an effective 
aperture of the time reversal mirror. We can rewrite the point spread function in terms of 
a normalized Gaussian as 

with the factor in front of the normalized Gaussian also equal to 



2ira e 

This means that when there is strong multipathing the integral of the point spread function 
over R d is not equal to one but to this ratio, which can be much smaller than one if a e a. 
Multipathing produces a tighter point spread function but there is also loss of energy, as of 
course we should expect. 

A more direct interpretation for the effective aperture can be given if the time reversal 
mirror has a Gaussian aperture function 

The point spread function W has now the form 



with 



a g = 



L 

ka„ 



and the effective aperture a g given by 



Clearly, a g w a e when there is strong multipathing and a e a. Written with a normalized 
Gaussian the point spread function for a Gaussian aperture has the form 



w(v) = - 



aj (27^/ 



2 
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Fig. 3.1. A directed field propagates from a distributed source of size a s toward the time reversal 
mirror of size a. The time-reversed, back-propagated field depends on the location of the mirror relative to 
the direction of the propagating beam. 



3.4. Broad-band time reversal for distributed sources. For a distributed source, 
its support a s is large compared to the F resn el number 6 so the ratio (3 = a_J6 is large. In 
this case we can compute the average of (3.7) the same way as we did in (3.8) and we find 
that 

(3.14) (* B (L,xo,£,i)> 

= (2n) d e^<- k ^M^/P) J (VF(L,x ,-^))e-* fecot g(- Co fc)^ 

ifpo-£-fcoc„t) ; rc/a\ f dydwc dk j(^m_ w . y _ fcc ^t) _dl 3 u 1 2 

= e HPo * *° Cot ^o(£//3) J (o^y+i *^(y)e ^ e 6 ff(-cofc)- 

The y integral on the right gives the Fourier transform of the aperture function xa{y) so 
with wo = cofco and a change of variable from fc to u> = cok we have 

(3.15) (* s (L,x ,e,t)) 



Here the star denotes convolution with respect to the spatial variables x, and a e is the 
effective aperture defined by ( [3.13; ). 

When multipathing is weak we can ignore the Gaussian factor in the convolution and 
we have 

(3.16) (* s (L,x ,£,i)) 

J 27T \UJ +LUJ 

In the opposite case, when there is strong multipathing and the effective aperture is much 
larger than the physical one, a e 3> a, we have 

(3.17) (v& B (L,x ,e,i)) 
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To interpret these results we note first that a distributed source function of the form 
( 3.11 ) can be considered as a phased array emitting an inhomogeneous plane wave, a beam, 
in the direction (k, pn), within the paraxial or parabolic approximation. The ratio |po|/fc 
is the tangent of the angle the direction vector makes with the z axis, and L|p |/fc is the 



transverse distance of the beam center to the center of the phased array (see Figure 3.1). If 



for each uj the beam displacement vector LcqPq / (ujo +u>) is inside the set A occupied by the 



time reversal array, then we recover at the source the full pulse in (3.16), time-reversed, 

(* s (L,x ,£,<)) = e^-^MWM-t)- 

If, however, for some frequencies the transverse displacement vector is outside the time 
reversal array, these frequencies will be nulled in the integration and a distorted time pulse 
will be received at the source. Depending on the position of the time reversal mirror relative 
to the beam, high or low frequencies may be nulled. 

In a strongly multipat hing medium the situation is quite different because the expression 



( 3.17 ), or more generally ( 3.15 ), holds now. Even if the beam from the phased array does not 
intercept the time reversal mirror at all, we will still get a time reversed signal at the source 
but with a much diminished amplitude. If the beam falls entirely within the time reversal 
mirror then the time reversed pulse will be a distorted form of g(—t), with its amplitude 
reduced by the factor (a/a e ) d . An interesting and important application of the time reversal 
of a beam in a random medium is the possibility of estimating the effective aperture a e 
by pointing the beam in different directions toward the time reversal mirror, measuring the 
time reversed signal that back propagates t o the source, that is, to the phased array, and 
inferring a e by fitting the measurements to ( 3.15| ). 



4. Summary and conclusions. We have analyzed and explained two important phe- 
nomena associated with time reversal in a random medium: 

• Super-resolution of the back-propagated signal due to multipathing 

• Self-averaging that gives a statistically stable refocusing 



Our analysis is based on a specific asymptotic limit (see Section 2.1) where the longitudinal 
distance of propagation is much larger than the size of the time-reversal mirror, which in 
turn is much larger than the correlation length of the medium, fluctuations in the index of 
refraction are weak, and the wave length is short compared to the correlation length. This 
asymptotic regime is more relevant to optical or infrared time reversal than it is to sonar or 
ultrasound. We have related the self-averaging properties of the back-propagated signal to 
those of functionals of the Wigner distribution. Self-averaging of these functionals implies 
the statistical stability of the time-reversed and back-propagated signal in the frequency 
domain, provided that the source function is not too broad compared to the Fresnel number 



(2.4). Time reversal refocusing of waves emitted from a distributed source is self-averaging 
only in the time domain. 

We apply our theoretical results about stochastic Wigner distributions to time reversal 
and discuss in detail super-resolution and statistical stability in section |^. 

Appendix A. The white noise limit and the parabolic approximation. 



We collect here some comments on the scaling analysis of section 2.1 and refer to [jl], |27J, 
[33| for additional comments and results on scaling and asymptotics in the high-frequency 
and white-noise regime. 



The dimensionless para mete rs 6, £,7 introduced by (2.2) in Section 2J, a long with the 



Fresnel number 9 defined by (2.4), lead to the scaled parabolic wave equation (2.5). If we do 
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not make the parabolic approximation and keep the ip zz term we have the scaled Helmholtz 
equation, with the phase e lkz removed, 



(A.l) 



2/)2 



e^0_ 

S 2 



■i/) gz + 2ik9ip z + 9 2 A x ip + 



k 2 6 



e 



0. 



Here, as in (|2.5| ), we relate the strength of t he fl uctuations a to s and 8 by (2.6). Is the 
parabolic approximation valid in the ordering (2/7), 6 -C £ <C 8 <C 1, that we have analyzed? 
The answer is yes but not before both 9 and e limits have been taken, in which case the 
scaled Wigner distribution (^^) conver ges to the Liouville-Ito process that is defined by the 
stochastic partial differential equation ( 2.19| ). 

It is in the white noise limit e — > 0, with Fresnel number 9 and 6 fixed, that the parabolic 
approximation is valid for (A.l), as was pointed out in Q. This is easily seen if the random 
fluctuations fi are differentiable in z. The parabolic approximation is clearly not valid in the 
high frequency limit 9 — > 0, before the white noise limit e — > is also taken. In the white 
noise limit, the wave function if>(z,x) satisfies an Ito-Schrodinger equation 



(A.2) 



2ik9d z ip + 9 2 A^dz 



ik 3 6 2 



Ro(Q)ipdz + k 2 6i>d z B(-, z) = 0. 



Here Rq is the integrated covariance of the fluctuations /1 given by ( pT5| ) and ( p,13[ ), and 
the Brownian field B(x,z) has covariance 

(B(x, zi)B(y, z 2 )) = i? (x - y)zi A z 2 . 



This Ito-Schrodinger equation is the result of the central limit theorem applied to ( A.l ). Let 



S e (x,z) 



1 



/z(x, -)ds. 

£ 



Then, as e — > this process converges weakly, under suitable h ypoth eses, to the Brown- 
ian field -B(x, z) with the above covariance. The extra term in ( |A.2| ) is the Stratonovich 
correction. 

The white noise limit for stochastic partial differential equations is analyzed in |l(J and a 
rigrous theory of the Ito-Schrodinger equation is given in 1 1 1 . The ergodic theory of the Ito- 
Schroedinger equation is explored in |17j . Wave propagation in the parabolic approximation 
with white-noise fluctuations is considered in detail in |l^, 30 j. 

The scaled Wigner distribution for the process ip, defined by ( |2.8| ), satisfies the stochastic 
transport equation 



(A.3) d z W e (z,x,p) 



V x W$(z,x, p)dz 



k 2 S 2 

Hp 

ikS 



c/q 



^-^i?o(q) I We{z,x,p+ — ) - W r e (z,x,p)J dz 



!M [JB- e i**/s (w e (zxp- 
29 J (2ir) d \vv eK z,x,p 



§)-»«,x,P + §) 



d z B(q, z), 



which is derived from ( A.2 ) equation using the Ito calculus. The Wigner process Wg con- 
verges in the limit 9 ^ to the Liouville-Ito process defined by the stochastic partial 
differential equation ( 2.1 9| ) . 
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Appendix B. Decorrelation of the Wigner process. 

B.l. Proof of Theorem \2.l[ We give here the proof of Theorems 2d and 2/2. We 
consider Theorem 2.1 first. It will follow from the Lebesgue dominated convergence theorem 
if we show that for pi ^ p 2 : 



(B.l) 



E {W s {z, x, Pl )W s {z, x, p 2 )} - E {W s (z, x, Pl )} E {W s (z, x, p 2 )} -» 



as 6 — > because the function Ws is uniformly bounded and E {Ws(z, x, pi)} does not 
depend on 5. Furthermore, the correlation function at the same spatial point but for two 
different values of the wave vector, Ug (z,x, pi,p 2 ) = E {Wg(z, x, pi)Ws(z, x, p 2 )} is the 



solution of ( p.lq ) with N = 2 and the initial data 

wj 2) (0,xi,pi,X2,p 2 ) = V^/(xi,pi)W A / (x 2 ,p 2 ), 

(2) 

evaluated at xi = X2 = x. Therefore Ul may be represented as 

c/f ( 2r ,x 1 , Pl ,x 2 ,p 2 ) = £;{w J (xJ(«),Pj(«))w>(x3( i f),pS(«))}. 

The processes X^' 2 and P^' 2 satisfy the system of SDE's (2.16) which may be more explicitly 
written as 

(B.2) 



dPl = - 



a(0)dB 1 (z) + -a 



X 5 " X 2 



dPj = - 



cr(0)dB 2 (z) + -a 



dX\ = -P}dz, 



dX 2 s = ~P 2 s dz 



dB 2 (z) 
dB 1 (z) 

with the initial cond ition s Xj' 2 (0) = x, P™(0) = p m , m = 1,2. H ere g 2 (0) = D, the 
diffusion coefficient ( 2.15 ) and the coupling matrix <r(x) is given by ( 2.17 ). Recall that 
Wg(z,x, p, k) — Ws(z, x, p/k; 1) and we need only consider the case k = 1. 

It is convenient to introduce the processes X 1 ' 2 and P 1 ' 2 that are solutions of ( 2.16| ) 
with no coupling: 



(B.3) 



dP m 

X l,2 



(0) 



-a(0)dB m (z), dX m = -P rn dz, 
= x, P m (0) = p m , m = l,2 



and define the deviations of the solutions of the coupled system of SDE's (B.2) from those 
of (p3h: Zf = Xf - X"\ Sf = Pf - P rn . Then we have 



(B.4) 



dSl =--a 



X 2 



dB 2 (z), dS 2 s 



1 

'2 a 



dB 1 (z) 



-Sgdz, 



dZ 2 



-S 2 s dz 



dZ\ 

with the initial data SJP(O) = Z m (0) = 0. Define 
(B.5) V(X 1 ,X 2 ,P 1 ,P 2 ,Zi,Z 2 ,Sj,S 2 ) 

= w^x 1 + zj, p 1 + s})Wi(x 2 

then we have with the above notation 
(B.6) 



Z 2 ,P 2 



S 2 5 )-W I (X 1 ,P 1 )W I (X 2 ,P 2 



E {W s (z, x, pi)W 5 (z, x, p 2 )} - E {Ws(z, x, Pl )} E {W 5 (z, x, p 2 )} 
= E {V{X\z) 1 X 2 (z) 1 P\z) 1 P 2 (z)Xs{z), Z 2 (z), Sj(*), S 2 (z))} 
<Ci?{|ZKz)| + |Z 2 (z)| + |Si(z)| + |S 2 (z)|} 
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since Wj is a Lipschitz function. 

Let us assume for simplicity that the correlation function -R(x) has compact support 
inside the set |x| < M. Then the coupling term in flB.2| ) is non-zero only when |X^_— X 2 | < 
MS. We introduce the processes Qs = P] — P 2 and Yg — Yl — y2 
satisfy the SDE's 



X^ that govern (B.4). They 



(B.7) 



dQs 



dB, dY s 



-Qsdz, 



Qi(0) = Pi-p 2 , Y 5 (0) = 



with B — B 1 — B 2 being a Brownian motion. 

In order to prove the theorem we show that the coupling term a(-) in (B.2) introduces 
only lower order correction terms, that is, S™ and Z™ are small. We show first that after a 
small 'time', t, the points X™ are driven apart since Qs(0) = Pj(0) — P 2 (0) 7^ 0. Then we 
show that after the points have separated the probability that they come close, so that the 
coupling term cr(-) becomes non-zero, is small. This "non-recurrence" condition requires that 
the spatial dimension d > 2. It follows that to leading order the points X™ are uncorrelated 
when d > 2 and that the coupling term introduces only lower order corrections. A similar 
argument for d = 1 would require an estimate on the time that points that are originally 
separated in the spatial variable spend near each other, where the coupling term in ( |B.2| ) is 
not zero. 

We need the following two Lemmas. The first one shows that particles that start at 
the same point x with different initial directions pi and P2, get separated with a large 
probability: 

Lemma B.l. Let Y$, Qs solve (B.7) with Ys{0) = 0, Q<j(0) = q ^ 0. Then for any 
e > there exists to(e) > that depends only on q = pi — P2 but not on S so that we have 



P 



|Yj(t)|>^) l - : Jo, all 7 ml.-) 



The second lemma shows that after the particles are separated, the probability that 
they come close to each other is small: 

Lemma B.2. Given any fixed r > and z > 0, ifYg, Qs solve \B.\ j with |Ya(0)| > r, 
Qs(0) =q^0, then P (M < s < z \Y 5 (s)\ < MS) -> as 5 -»■ 0. 

We prove Theorem 2.1 before proving Lemmas B.l and B.2: 

Proof. Let z and q = pi — P2 be fix ed an d defined as a bove. Given e > 0, then for any 
t < To(e) (with ro as defined in Lemma B.l), Lemma B.2 and the Markov property of the 
Brownian motion imply that 



/'(Sr':) = S:; ( r) |Y tf (r)|>^) I : 



and 



r(Z,rc) = Z-(r) + (:- > \Y 6 (r)\ > T -f ) I 



for (5 < So(r,e). Furthermore, 



(B.8) 



E{ \Z\{r)\ + |Z 2 (r)| + |Sj(r)| + |S 2 (r)| |Y a (r)| > 



< E{\Z](r)\ + |Z 2 (r)| + \S](r)\ + |S 2 (r)|} /(l - s) < Cr 
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because the function <r is uniformly bounded. Therefore we have 
£{V(X\X 2 ,P\P 2 ,Zi,Z2,Si,S2)} 

(B.9) =^|v(X 1 ,X 2 ,P 1 ,P 2 ,Zi,Z|,Sj,S|)||Y 5 (r)| > (\Y s {r)\ > ^ 

+ J e|v(X 1 , X 2 , P\ P 2 , Zj, Z 2 , Si, S 2 )| |Y,(r)| < ^} P (|Y,(r)| <^j=I + II. 

The second term above is small because the probability for Y$ (t) to be very small is bounded 
by Lemma B.l. More precisely, given e > and r < tq(e), Lemma |B.l| implies that 



(B.10) 



i7 < Ce. 



The first term in (B.9) corresponds to the more likely scenario that Yg at time t has 
left the ball of radius r|q|/2. We estimate it as follows. The probability that Y,5 re-enters 
the ball of radius M5 is small according to Lemma B.2. Moreover, if Yg stays outside this 
ball, the difference variables Z m and S m are bounded in terms of their values at time r. The 
latter are small if r is small. More precisely, using (B.8) we choose r so small that 

r\q\ 



EUZ s (t) 



|Z 2 (r)| + |Si(r)| + |S 2 (r)| |Y ff (r)| > 



< £. 



Then we obtain 



I < E jv(X\X 2 , P\ P 2 , Z], Z 2 , Si, S 2 )||Y 5 (r)| > iM J 

<E\v(X\X 2 ,P\P 2 ,Zl,ZlSl,Sj) \Y s (t)\ > ^ and inf |Y 5 ( S )| < Ms) 



brfjY,( a )|<M*||Y 5 (r)|>l^ 



|Z 2 (z)| + |Sj(*)| + |S 2 (z)| |Y,(r)| > ^ and jpf | Ytf ( a )| > Ms) 

Z T<S<Z J 



xP[ inf \Ys(s)\>MS\Y s (r)\>^ 

\ r<s<z Z 



= h + h. 



The term Ii goes to zero as S — > by Lemma B.2. However, if the conditions in / 2 hold, 
then 

Sf(z) = Sf(r), Zf(z) = ZT(r) - \{z - r)Sf (r). 



Therefore the term I2 may be bounded with the help of ( |B.8| ) by 
I % < E{\Zl(z)\ + |Z 2 (z)| + \S](z)\ + |S 2 (z)|||Y,(r)| > ?M and ^ \Y 5 (s)\ > M<j} 
< C7£;||Zj(r)| + |Z 2 (r)| + |Sj(r)| + |S 2 (r)| | |Y 5 (r)| > lM| < Cr. 

Putting together ( |B.9| ), ( p3 ■ 10| ) and the above bounds on I\ and I2, we obtain 
E{\ZUz)\ + \ZUz)\ + \SUz)\ + \Sj(z)\}<Ce 



for (5 < J and Theorem follows from (B.6). □ 
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B.2. Proof of Lemmas B.l and B.2 . We first prove Lemma B.l 
Proof. We write 



Qs(z) = q 



<r(0) - ^r(Y(s)/5) ) ,/Bi.s) = q+Q,,Ci 



Y s (t) = qt / Q s (s)ds 



so that 

Then we have 
(B.ll) 

and hence 



P (|Y 4 (r) + rq| > rr) < P sup |Q,(a)| >r) <Ct/ 

\0<S<T 

We let r = |q|/2 in the above formula and obtain 

rlq 



P( sup |Q s (s)| > r) <Cr/r 2 

\Q<S<T 



P |Y 4 (r)| < 



C 

|q| 2 



and the conclusion of Lemma 



B.l follows. □ 



Finally, we prove Lemma B.2 



Proof. Let rg be the first time Yg[z) enters the ball of radius MS: 

r A -=inf{z: \Y 5 (z)\<MS}, 

with Y s (0) = YV 0. For < a < 1 let Az = S 1 "", n = \z/Az], Ji = (iAz, (i+ l)Az) and 
p < 1. Note that until the time r 4 the process (Yg, Qg) coincides with the process (Y, Q) 
governed by (|B.7|) without the coupling term a(Yg/S). We find 



P(ts <z) <^ jp (I Y(iAz)| <M<S P ) + P^inf |Y(s)| < |Y(tA«)| > MS p ^j | . 

The process Y(s) is Gaussian with mean Y° and variance 0(s 2 ). Therefore, there is a 5 > 
such that for S < 5 



If we assume 
(B.12) 
then also 



P(|Y(iAz)| < MS P ) < CS dp . 
p < 1 — a 



P(t s <z)<nC (s dp + p( sup | Y(s) - Y° | > M[S P - S\]] 

V \0<s<Az // 

<C{5 dp+a - x + 5 a - x p( sup |B(s)| > M[S p -S]/Az)) 

\0<s<Az 

ja _ 1 g{B(Az) 2r }Az 2 '- 



< C 



(<Jp - £) 2r 

+a — 1 ja— 1— rp+3T*(l— a)/2 
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Note that with p < 1 — a and r large enough, there is a q > so that 

P(t S <z)< C8 q 



if d > 2 and Lemma B.2 follows. □ 

B.3. Proof of Theorem |2.2| . We need to show first that 

(B.13) Mz,*) = J Ws(z,x,p)dp 

is finite with probability one. The stochastic flow (X$(t,x, p),P$(t, x, p) is continuous in 
(£,x, p) with probability one, so Wg(z,x,p) — Wi(X.s(t, x, p), P$(f, x, p)) is bounded and 
continuous. It is, moreover, non-negative if Wi > 0. We know that 

E{W s (z, X ,p)}dp 

is finite and independent of <5, and the order of integration and expectation can be inter- 
changed by Tonelli's theorem. This theorem implies in addition that Jg(z,x.) is finite with 
probability one. 

We can now consider 

E{J$(z,x)} = J E{W s (z,x,p 1 )W s (z,x,p 2 )}dpidp 2 . 
The integrand is bounded by an integrable function uniformly in S because 

E{W s (z,x, Pl )W s (z,x,p 2 )} < E l ^{Wi(z,x,p 1 )}E 1 ^{W s {z,x,p 2 )h 

the right side does not depend on <5, and is integrable. Therefore by the Lebesgue dominated 
convergence theorem and the results of the previous Section we have that 

lim E{jUz,x)}=E 2 {J s (z,x)} 
and the right side does not depend on S. This completes the proof of Theorem |2.2|. 
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